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We investigate the critical properties of the spin-3/2 Blume-Capel model in two dimensions on a 
random lattice with quenched connectivity disorder. The disordered system is simulated by applying 
the cluster hybrid Monte Carlo update algorithm and re-weighting techniques. We calculate the 
critical temperature as well as the critical point exponents 7/1^, l3/iy, ajv, and v. We find that, 
contrary of what happens to the spin-1/2 case, this random system does not belong to the same 
universality class as the regular two-dimensional ferromagnetic model. 
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INTRODUCTION 

Experimental studies of the critical behavior of real 
materials are often confronted with the influence of im- 
purities and inhomogeneities |l|. For a proper interpre- 
tation of the measurements it is, therefore, important to 
develop a firm theoretical understanding of the effect of 
such random perturbations. In many situations the typi- 
cal time scale of the thermal fluctuations in the idealized 
"pure" systems is clearly separated from the time scale 
of the impurity dynamics, such that to a very good ap- 
proximation the impurities can be treated as quenched. 
The importance of the effect of quenched random disor- 
der on the critical behavior of a physical system can be 
classified by the specific heat exponent of the pure sys- 
tem, a^uTf^- The criterion due to Harris [3 asserts that for 
ttpttre > quenched random disorder is a relevant pertur- 
bation, leading to a different critical behavior than in the 
pure case (vifhich is the case of the three-dimensional Ising 
model). In particular, one expects [3| in the disordered 
system that v > 2/D, where i' is the correlation length 
exponent and D is the dimension of the system. Assum- 
ing hyper-scaling to be valid, this implies a = 2 — 0^ < 0. 
On the other hand for apure < disorder is irrelevant (as 
is the case of the three-dimensional Heisenberg model) 
and, in the marginal case apure = 0, no prediction can 
be made. For the case of (non-critical) first-order phase 
transitions it is known that the influence of quenched 
random disorder can lead to a softening of the transition 
y . Recently, the predicted softening effect at first-order 
phase transitions has been confirmed for 3D chstate Potts 
models with q > 3 using Monte Carlo |a, la H and high 
temperature series expansion [g] techniques. The over- 
all picture is even better in two dimensions (2_D) where 
several models with am.rp. > 0, Efl llll UJj and the 
marginal {apure = 0) [ij, |lj, |la, lla ll^ have been inves- 
tigated. 

In this paper we study another type of quenched ran- 
dom disorder, namely connectivity disorder, a generic 



property of random lattices whose local coordination 
number varies randomly from site to site. Specifically, 
we consider 2D Poissonian random lattices of Voronoi- 
Delaunay type, and performed an extensive computer 
simulation study of a Blume-Capel model. We concen- 
trated on the close vicinity of the transition point and 
applied finite-size scaling (FSS) techniques to extract the 
exponents and the "renormalized charges" f/2 and C/|. 
To achieve the desired accuracy of the data in reasonable 
computer time we applied the single-cluster hybrid algo- 
rithm ISJ to update the spins and furthermore made ex- 
tensively use of the re- weighting technique [19|] . Previous 
studies of connectivity disorder focusing mainly on 2D 
lattices have been realized by Monte Carlo simulations 
of g-state Potts models on quenched random lattices of 
Voronoi-Delaunay type for a = 2 20. 21.i2i|. g = 3 |2^ 
and q = 8 |24 L l25 l . In particular, it has been shown that 
for 9 = 2 jia |3i 113 and g = 3 nil the critical expo- 
nents are the same as those for the model on a regular 2D 
lattice. This is indeed a surprising result since the rele- 
vance criterion of the Delaunay triangulations reduces to 
the well known Harris criterion such that disorder of this 
type should be relevant for any model with positive spe- 
cific heat exponent 26]. This means that for q — 3, where 
apure > 0, one would expect a different universality class. 
On the other hand, for the present spin-3/2 model, where 
apure = 0, we show that the exponents indeed change in 
the Voronoi-Delaunay lattice type, turning out the situ- 
ation still more bizarre . In the next section we present 
the model and the simulation background. The results 
and conclusions are discussed in the last section. 



MODEL AND SIMULATION 

The Voronoi construction or tessellation for a given set 
of points in the plane is defined as follows ^27] . Initially, 
for each point one determines the polygonal cell consist- 
ing of the region of space nearer to that point than any 
other point. Then one considers that the two cells are 



neighboring when they possess an extremity in common. 
From the Voronoi tessellation the dual lattice can be ob- 
tained by the following procedure: (a) when two cells are 
neighbors, a link is placed between the two points located 
in the cells; (6) From the links one obtains the triangula- 
tion of space that is called the Delaunay lattice; (c) The 
Delaunay lattice is dual to the Voronoi tessellation in the 
sense that points corresponding to cells link to edges, and 
triangles to the vertices of the Voronoi tessellation. 

We consider now the two-dimensional spin-3/2 Blume- 
Capel model on this Poissonian random lattice. The 
Blume-Capel Model is a generalization of the standard 
Ising model |2g and was originally proposed for spin-1 
to account for first-order phase transition in magnetic 
systems |22l l2i| ■ The Hamiltonian can be written as 



C{K) = K^N[< e^ >-<e >% 



(3) 



H 



-jJ2 s,s, + Aj2s:, 



(1) 



<ij> 



where the first sum runs over all nearest-neighbor pairs 
of sites (points in the Voronoi construction) and the 
spin-3/2 variables Si assume values ±3/2, ±1/2. In eq. 
(1) J is the exchange coupling and A is the single ion 
anisotropy parameter. The second sum is taken over the 
N spins on a D-dimensional lattice. The case where S = 
1 has been extensively studied by several approximate 
techniques in two- and three-dimensions and its phase 
diagram is well established [H El HI El Hi El ■ 
The case S > 1 has also been investigated according to 
several procedures El El El El E El 113 ■ 

The simulations have been performed for A = 0, which 
is the simplest case, on different lattice sizes comprising a 
number N = 1000,2000,4000,8000, 16000 and 32000 of 
sites. For simplicity, the length of the system is defined 
here in terms of the size of a regular lattice L — N^^^ . For 
each system size quenched averages over the connectiv- 
ity disorder are approximated by averaging over i? = 100 
(N = 1000 to 4000), R = 50 (N ^ 8000) and i? = 25 
{N — 16000 and 32000) independent realizations. For 
each simulation we have started with a uniform config- 
uration of spins (the results are however independent of 
the initial configuration). We ran 2.52 x 10^ Monte Carlo 
steps (MCS) per spin with 1.2 x 10^ configurations dis- 
carded for thermalization using the "perfect" random- 
number generator [43| . We have employed the hybrid 
algorithm J18j where we included n Wolff clusters (here 
n = 5) intercalated by one Metropolis single-spin flip 
sweep. This algorithm has been shown to be quite effec- 
tive for spin-3/2 models [ll|. For every 12th MCS, the 
energy per spin, e — E/N , and magnetization per spin, 
m — ^j Si/N , were measure and recorded in a time se- 
ries file. 

From the series of the energy measurements we can 
compute, by re- weighting over a controllable temperature 
interval AT, the average energy and specific heat 



where K — J/ksT, with J — 1, and ks is the Boltz- 
mann constant. In the above equations < ... > stands 
for thermodynamic averages and [..^av for averages over 
the different realizations. Similarly, we can derive from 
the magnetization measurements the average magnetiza- 
tion, the susceptibility, and the magnetic cumulants. 



m{K) = [< \m\ >]av, 
X{K) = KN[< m^>- < \m\ >\ 



U2{K) = [1 



Ua{K) = [1 



<m^ > 



3 < \m\ >'^ 



<mr > 



3 < \m\ >^ 



(4) 



(5) 



(6) 



(7) 



Further useful quantities involving both the energy and 
magnetization are their derivatives 



d[< \m\ >] 
dK 



[< \m\E > - < \m\ >< E >]av, (8) 



dln[< \m\ >]av _ r< \m\E > 



dK 



dln[< |m^| >]fl 
dK 



< TO > 



,< \m'^\E > 



< 771^ > 



-<£;>! 



(9) 



< E >U. (10) 



In the infinite-volume limit these quantities exhibit sin- 
gularities at the transition point. In finite systems the 
singularities are smeared out and scale in the critical re- 
gion according to 



C^Creg+L-''^''fc{x)[l + ...], 



(11) 



[<\m\>U^L-^/''f,nix)[l + ...], (12) 

X = L-^/Vx(2^)[l + ■■■], (13) 



dln[< \m\P >]a.y 
dK 



- rV^ 



fp{x)[l 



(14) 



u{K) ^[<E >U/N, 



(2) 



where Creg is a regular background term, i/, a, /3, and 
7 are the usual critical exponents, and fi{x) are FSS 
functions with x — {K — Kc)L^^'^ being the scaling vari- 
able, and the brackets [1 -I- ...] indicate corrections-to- 
scaling terms. We calculated the error bars from the 
fluctuations among the different realizations. Note that 
these errors contain both, the average thermodynamic 
error for a given realization and the theoretical variance 
for infinitely accurate thermodynamic averages which are 
caused by the variation of the quenched, random geome- 
try of the lattices. 



RESULTS AND CONCLUSION 



regular lattice exponent v = I. 



By applying standard re-weighting techniques to each 
of the R time-series data we first determined the temper- 
ature dependence of Ci{K), Xi{K),..., i = 1,...,R, in the 
neighborhood of the simulation point Kq. Once the tem- 
perature dependence is known for each realization, we 
can easily compute the disorder average, e.g., C{K) — 
X]i=i Ci{K)l Ri ^^^d then determine the maxima of the 
averaged quantities, e.g., CmaxiKmax) = maxKC'{K). 
The variable R represents the number of replicas in our 
simulations. 

In order to estimate the critical temperature we calcu- 
late the second and fourth-order Binder cumulants given 
by eqs. (7) and (8), respectively. It is well known that 
these quantities are independent of the system size and 
should intercept at the critical temperature [4^. In 
Fig. ^ the fourth-order Binder cumulant is shown as 
a function of the K for several values of N. Taking 
the largest lattices we have K^, — 0.1844(1). To esti- 
mate C/4 we note that it varies little at K^, so we have 
C/| = 0.482(6). From the second-order cumulant we simi- 
larly get Kc = 0.1845(1) and [/| = 0.579(8). One can see 
that the agreement of the critical temperature is quite 
good and U^ is definitely far from the universal value 
f/| ~ 0.61 for the same model on the regular 2D lattice. 
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FIG. 1: Fourth-order Binder cumulant as a func- 

tion of K for several values of the system size N = 
1000, 2000, 4000, 8000, 16000 and 32000. 

The correlation length exponent can be estimated from 
the derivatives given by eq. (15). Figure [21 shows the 
maxima of the logarithm derivatives as a function of the 
logarithm of the lattice size L ior p — 1 and p = 2. 
From the linear fitting one gets v = 0.85(2) {p = 1) and 
z^ — 0.917(8) { p — 2), which is again different from the 



7.8 



:7.2 






^ 6.8 



V 

:=, 6.2 






5.2 




3.2 



3.8 



4.2 

LnL 



4.8 



5.2 



FIG. 2: Log-log plot of the maxima of the logarithmic deriva- 
versus the lattice size L — N^'^ for p — 1 
2 (square) . The solid lines are the best linear 



tive 



dln[<|m|P>] 



dK 

(circle) and p ■- 
fits. 



In order to go further in our analysis we also com- 
puted the modulus of the magnetization at the inflection 
point and the maximum of the magnetic susceptibility. 
The logarithm of these quantities as a function of the 
logarithm of L are presented in Figures O and 01 respec- 
tively. A hnear fit of these data gives /3/i^ = 0.331(9) 
from the magnetization and j/v = 1.467(9) from the 
susceptibility which should be compared to /3/j^ = 0.125 
and j/i' — 1.75 obtained for a regular 2D lattice. 
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FIG. 3: Plot of the logarithm of the modulus of the magne- 
tization at the inflection point as a function of the logarithm 
of L = N^'"^ . The solid line is the best linear fit. 
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FIG. 4: Log-log plot of the susceptibility maxima Xmax as a 
function of the logarithm of L = N^'^ . The solid line is the 
best linear fit. 



The specific heat can also be analysed in this case but, 
as it happens in other models |2lJ,|23j, we cannot find a 
clear unambiguous support for a definite scaling. Figure 
|S1 shows the maximum of the specific heat Cmax as a 
function of L. Least-squares fits to a logarithmic Ansatz 
Cma. ^ Bo + BilnL give Bo = 0.44(6), Bi = 0.72(1) 
and is shown by the full line in figure |31 The dashed line 
in this figure corresponds to a pure power-law Ansatz, 
C„„^ = cL"'/" with c = 1.475(5) and a/u = 0.202(5). 
From these results one can slightly see a better agreement 
with the logarithmic Ansatz. 
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FIG. 5: Specific heat maxima Cmax as function of L = N^''^ . 
The solid line is the best fit to an a ~ (Log) Ansatz and 
the dashed line to a power law Ansatz. 



Thus, from the above results, there is a strong indica- 
tion that the spin-3/2 Blume-Capel model on a Voronai 



lattice is in a different universality class than its regu- 
lar lattice counterpart. This poses, in addition to the 
q — i Potts model in two dimensions, and taking into ac- 
count the extensive study done Janke and Weigel on the 
Harris-Luck criterion for random lattices |26|, another 
open question to be answered in more general terms. 
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